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Developments in QCD at finite density are reviewed. I begin by discussing some new algorithms 
which have been applied to other theories with sign problems. Then I discuss the method of 
analytic continuation in QCD using a series expansion and review some of the results obtained 
using this method. By now there are several different simulations using the method, and together 
they give estimates of the systematic lattice effects, which turn out to be controlled. Finally 
I discuss a direct comparison of some of these lattice predictions with new experimental data 
which results in a very pleasant agreement. 



The XXVIII International Symposium on Lattice Field Theory, LatticelOlO 
June 14-19, 2010 
Villasimius, Italy 



* Speaker. 



© Copyright owned by the author(s) under the terms of the Creative Commons Attribution-NonCommercial-ShareAlike Licence. 



http://pos.sissa.it/ 



QCD with fiy^Q 



Sourendu Gupta 




Figure 1: A conjectured phase diagram of QCD for Nf = 2 and for Nf = 2 + 1 when the strange quark mass 
is not much smaller than Aqcd- In the chiral limit there is a tri-critical point, from which emerges an Ising 
critical line whose intersection with the plane of physical quark mass is the QCD critical point. 

1. Introduction 

QCD at finite baryon density is interesting because of two reasons: first that tliere is a program 
of experimental studies covering five colliders, running and planned, wliich will look at this prob- 
lem, and second, that it does not seem open to standard methods of attack in lattice gauge theory 
due to a sign problem [|l|]. In this review I will bring together evidence that the problem is still 
open to a fruitful attack using small modifications of the usual tools of lattice gauge theory, and 
give some of the main physics results. The context of these first results is the conjectured phase 
diagram of Figure |l| [^. 

Any Monte Carlo integration process suffers from a sign problem if the integrand is not real 
and positive definite. For the QCD action with a chemical potential on the baryon number, the 
determinant of the Dirac operator, which is the quark part of the measure, obeys the condition 

det(D + m + ;U7o)* =det(D + m-^*7o), (1.1) 

where D is the massless Dirac operator, m is the mass, jx is the baryon chemical potential, and * 
denotes complex conjugation. For any generic complex chemical potential this shows that there 
is a sign problem. For pure imaginary jj. (including n = 0), the determinant is real, and one can 
further prove its positivity by considering its commutation with 75. 

This sign problem is not necessarily mild. Baryonless random matrix theory seems to predict 
that for /I < the distribution of signs is Gaussian and becomes Lorentzian at larger [|3|]. In 
either case the problem is severe. An earlier work had estimated the contours of the variance of 
the phase of the quark determinant and found that this decreases at high temperatures, where the 
problem could therefore become easier 
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Figure 2: The line is the phase boundary of the 3d XY model, found using the worm algorithm [||]. The 
instability of the complex Langevin method is illustrated in this figure [pi]]: this method converges to the 
wrong probability measure. The difference between the actions for the true and converged distributions is 
colour coded. 

This review is structured as follows. Section ^ presents an overview of very interesting new 
attempts to attack the sign problem directly; unfortunately they are not yet at the stage where 
they can be applied to QCD. Section |3| reviews the Madhava-Maclaurin series expansion^ method 
which has yielded first results on the phase diagram and on some other measurable quantities. 
Finally in Section ^ first results from experiments are reported along with comparisons to lattice 
QCD predictions. 

2. Trial algorithms 

The class of algorithms which has had the most attention till now is reweighting: perform the 
Monte Carlo procedure at a point in the phase diagram where there is no sign problem, and then find 
expectation values of operators by choosing an appropriate weight for each configuration. Various 
problems with this process are by now well-known; they include large errors due to cancellations 
and inaccurate sampling, which become exponentially large with the volume. The Budapest group 
applied this method to the problem of locating the QCD critical point |^. The Bielefeld-Swansea 
algorithm is a variant of this method which expands the determinant in a series in jj. [^. There 
have been no developments in this class of algorithms since it was reviewed in 2008 [^]. 

Two new classes of algorithms are being tested currently, and, although their applicability to 
QCD is not yet clear, they are interesting enough to merit some discussion. Interestingly, they 
become easy to compare because the two algorithms have been, deliberately, applied to the same 
model recently. This is the 3-d XY model at finite chemical potential, which has the action 

S = -j8 ^cos (0, - 0,+^ - /At5^,f) . (2.1) 

'in the 14th century Madhava of Sangamagrama developed the series expansion for functions and estimates of the 
error terms which later came to be associated with the name of Maclaurin. 
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This suffers from a sign problem when /i 7^ 0. 

One approach exploits the fact that sign problems are not inherent to the physics of a 
system, but to specific representations. By a clever transformation of fields which amounts to 
redefining the theory in terms of fluxes of particles along links, they reduce it to a form without 
a sign problem, although the theory then becomes non-Iocal. However, in this representation it 
becomes amenable to a numerical attack using the "worm algorithm" This work then sets out 
a finite-size scaling theory which describes the point at which it becomes energetically favourable 
to add one more particle to the ground state. The simulation results allow the extraction of finite 
size scaling parameters which can then be used to determine the phase diagram. 

The other approach resurrects an old idea — the complex Langevin method, wherein one ad- 
dresses sign problems by complexifying the fields while the noise remains real. Earlier works 
had been plagued by runaway directions and associated numerical instabilities, now brought under 
control by the use of adaptive step-size integrators. For a while a proof of convergence of such 
methods seemed to be within reach [|l0|]. However, it turns out that there may be a convergence to 
the wrong result [11]. This is illustrated in Figure ^, which shows that the problem arises mainly at 



small temperature and large chemical potential. Since this region is similar to that in which QCD 
has large sign fluctuations [^], a better understanding of the origin of this problem may throw light 
on applications to QCD. 



In the next section we turn to the algorithm, first described in [13], which is now used by many 



groups, and has begun to yield many consistency tests and, possibly, even contact with experiment. 

3. The Madhava-Maclaurin series expansion 

The pressure of QCD matter in a grand canonical ensemble can be expanded in a Madhava- 
Maclaurin series around the point = to obtain 

P(r, At) = PiT) + ^X^'\t) + ^X^'Kt) + ... (3.1) 

where all the coefficients are computed at = 0. P{T) is the pressure at zero chemical potential, 
X^^\T) is called the quark number susceptibility (QNS) [12] and all the X^'^\T) are generically 



called non-linear susceptibilities (NLS). It was suggested that the NLS could be measured in = 



simulations, and the feasibility was demonstrated by computations in quenched QCD 013|]. More 
recently, within the last year, there have been attempts to compute these coefficients by simulating 
QCD at imaginary chemical potential and fitting extrapolating functions to the data [Q] (we will 
return to a discussion of this later). 

3.1 Computational effort 

The ;f^"^'s are combinations of quarks loops with insertions of 70 up to n times [15]. These 



quark loop traces are obtained through stochastic noise averages. One measure of the feasibility of 
such measurements is to examine the signal to noise ratio in the measurements when the number 
of noise vectors is Ny, i.e., the ratio of the mean and square root of the variance of such a trace in 
one configuration. When the ratio is large, the measurement is easy. Such a measure was reported 
using staggered quarks on a 4 x 24^ lattice at T /Tc = 0.75 and Nf = 2 when the quark mass is 
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Figure 3: Signal to noise ratios for various fermion traces which enter the evaluation of the NLS (see [Op 
for explanations of the notation). Circles denote data for staggered quarks Jl^ and boxes for P4 quarks []16|]. 
The red circles denote measurements of quantities which should be exactly zero. 



tuned to give nij^ = 230 MeV [15]. Here we add results from P4 action with Nf = 2 + \, with light 



quark masses tuned to the same value of niji and r/r^. = 0.84 [16]. In both cases Ny = 400 and the 



signal to noise ratios are comparable (see Figure |3|). A direct comparison with Asqtad quarks is not 



available, but from the claim that 50% of the noise in pi 70 is due to stochastic estimators, one finds 
that the signal to noise ratio for that Dirac operator is comparable. 

Near Tc autocorrelations between successive configurations is large — of the order of 200-250 
trajectories. Assuming that it takes about 200 fermion matrix inversions per trajectory, and that we 
use A^v = 500 for every decorrelated configuration, then, since it takes 18 inversions per measure- 



ment (up to the 8th order of the expansion in eq. p.l\ ), the ratio of CPU times for a measurement 
to that for generating a decorrelated configuration is 0.24. The marginal cost of measurement is 
small. Well inside the hot phase, at 27^, the autocorrelation time drops to about 4 trajectories, 
whereas A^v = 100. The ratio of CPU times for measurement to generating decorrelated configu- 
rations climbs to 4.5, however, with relatively small expenditure of CPU time. As a result, direct 
measurements of the NLS are highly feasible. An added attraction is that configurations which 
have been generated for any finite temperature study can be reused for such analysis, thus reducing 
the marginal cost even further. 

3.2 Series Analysis for the critical point 

Methods for analysis of series expansions of the free energy or its derivatives are well-known 
in statistical mechanics, and have been used successfully in many cases [[l^]. While the core of the 
analysis is the same, there are interesting differences between these older works and the application 
to QCD, which lead to differences in the method of analysis [ p3| , [l9| ]. The series coefficients in the 
older works came from exact enumeration of graphs, corresponding to infinite volume systems, so 
that a coefficient was either known exactly or not known. In the present case the series coefficient 
is evaluated on finite lattices with statistical errors through a Monte Carlo process. 
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Figure 4: Estimates of the radius of convergence of the series expansion in eq. from finite temperature 
simulations at 7" = 0.94Tc with two flavours of staggered quarks and the bare quark mass tuned to give 
m^: = 230 MeV and lattice spacing 1 / {6T). Open symbols correspond to estimates of z* and filled symbols 
to T„- When LT increases from 2 to 4 «*(L) increases significantly, as shown. 



The main point of the analysis is that the series for the quark number susceptibility can be 
analyzed for its radius of convergence. The series for the QNS is 

2!!^ = ^^;(I)4,(«(r)4r=,(«,r) + ... (3.2, 

where z = }i/T and each of the dimensionless combinations T'^^^x^'^\T) is the direct output of a 
lattice computation. Since the expansion is in z at fixed value of T, it is equivalent to an expansion 
in jU. If the divergence occurs at a temperature at which all the series coefficients are positive, then 
the non-analyticity occurs for real values of z, and the divergence can be identified with the critical 
point of QCD. 

The radius of convergence can be found by several methods, all of which correspond to com- 
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paring the series against another with a known singularity. The best known definitions are 

l/(«-2) 



and z„ = — T^rT;;^ ■ (3-3) 



(« + 2)!r"-4;t:W' l«!;^{2)/r2 



The star and bar do not indicate complex conjugation. The common limit as « — > oo of both is the 
radius of convergence of the series. This test is closely coupled to a finite volume scaling analysis. 

The reason is the following. If there is a critical point at some (/x^ , T^) then the QNS diverges 
there on an infinite volume system. However, on any finite volume, L^, there may be a peak, but 
no divergence. As the system size decreases, the peak becomes broader and lower. As a result, z* 
and z* may seem to give a finite radius of convergence for n < n^{L). For larger n both z* and z* 
will then become larger and larger, since there is no actual divergence in the series for the QNS. 
As L increases, one would find (L) also increasing without limit. In simulations of QCD with 2 
flavours of staggered quarks with the bare quark mass tuned to give m-^ = 230 MeV, such behaviour 
is actually seen at one temperature (see Figure ^. At this temperature all the NLS are positive, so 
the limiting singularity is at real values of pL. We can then identify such a temperature with and 
the corresponding estimate of the radius of convergence with pL^ /T^ . 

Such estimates of the critical point have been made with two different lattice spacings using 
two flavours of staggered quarks on large volume lattices A computation with 2+1 flavours 

of P4 quarks at almost the same value of rriji has also been performed with large lattices [ |20| ] and 



preliminary estimates of the radius of convergence have been reported [|21[]. These are collected 
in Figure ^. Since large volumes are crucial to obtaining a stable estimate of the critical point, 
older computations with smaller volumes have not been added into this figure even if they have 
realistic values of niji. A by-product of this choice is that all the points in the figure use the same 
computational technique, albeit with different lattice spacings and quark actions. 

Interest in the critical point is enhanced because of the possibility that heavy-ion collision 
experiments may observe it. Fireballs produced in such collisions undergo chemical freezeout at 
certain values of jx and T which change with the center-of-mass energy, \/S, of the colliding nuclei. 
The chemical freezeout point is relevant if one tries to use fluctuations of conserved quantum 
numbers as probes of the critical point: we shall return to this argument later. The freezeout curve 



is parametrized in [|24|]. This has been superposed on the phase diagram in Figure g by using the 
scale Tc = 175. 

A pleasant fact emerges from Figure ||: that lattice spacing effects can be bounded in magni- 
tude by currently available computations. The difference between different kinds of actions is, of 
course, a finite lattice spacing artifact. The magnitude of the lattice spacing effect estimated from 
two different spacings with the same action turns out to be comparable with that from a comparison 
of two different actions at nearly the same lattice spacing. 

Interestingly, the effect of the strange quark on the end point seems to be under control. It has 
long been known that in the Columbia plot, the physical point corresponds to a thermal cross over 



[22], as a result of which the topology of the phase diagram of realistic 2-1-1 flavour QCD is the 
same as for two flavours, as in Figure [T[ As the strange quark mass is increased and the light quark 
mass is reduced, the thermal crossover passes through a critical point into a first order transition. 
It turns out that the line of critical points is far from the physical point: the pion needs to be about 
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1.1 




Figure 5: Estimates of the critical point of QCD from lattice computations with niji ~ 230 MeV and Lnijc > 
4. The points in pink are obtained from computations with 2 flavours of staggered quarks and the point in 
brown from 2+1 flavour P4 quarks. The values of ji and T along the freezeout curve are parametrized by 
and have been used in conjunction with Tc = 175. Some values of VS have been marked along the 
freezeout curve. 
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Figure 6: The variation of the radius of convergence with in partially quenched computations with the 
staggered sea quark mass tuned to give m„ ~ 230 MeV. As the valence quark mass is changed the partially 
quenched pion mass is m^'^. 
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Figure 7: Using the series expansion to describe the data obtained through direct simulations at imaginary 
chemical potential. 



10 times lighter and the strange quark mass about 3 times lighter [g3|]. This could be one reason to 
suspect that the numerical values for the critical end point in 2 flavour and 2+1 flavour QCD may 
not be very different. Such an argument is compatible with the results collected in Figure ||. 

The major remaining effect is due to the light quark mass being larger than physical. The only 
exploration of this effect till now is a partially quenched computations with a = 1/{4T) and two 
flavours of staggered quarks with a sea quark mass tuned to give Mj^ = 230 MeV [|25|]. Interestingly, 
an interpolation between the measured value of the radius of convergence was consistent with the 
result of an unquenched computation with P4 quarks tuned to give ruji = 550 MeV [Q]. In Figure 
|6|the earlier results are extended by adding a similar analysis for a = 1/ {6T). The extrapolation to 
the physical value of niji shows a 15% drop in the value of jj.^. 



3.3 Extrapolation of observables 

Apart from the prediction of the critical point, the series expansion could be used to extrapo- 
late measurements to finite chemical potential. Tests of such extrapolations are whether they can 
describe measurements made directly through simulations at imaginary chemical potential. The 



most straightforward extrapolation is to use the series. A preliminary attempt [27] is shown in Fig- 
ure ^ One sees that adding new terms in the series improves the extrapolation only marginally in 
jj.. Closely related to this exercise are attempts to extract the series coefficients from measurements 
obtained in direct simulations at imaginary chemical potential. It was shown recently JT4| ] that 
simple series descriptions of the data obtained at finite imaginary chemical potential are inefficient 
and more complicated forms are needed to perform the extraction of the series coefficients. 

Clearly, when the series expansion of a physical quantity is close to divergence, then a trun- 
cated sum is clearly not the best way to find the value of this quantity at finite jj.. One must search 
for series resummation techniques. One method which has been widely used for resummation of 



high temperature series of spin models [|18[] is to determine Pade approximants using them. There 
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Figure 8: Extrapolation of the measurement of the QNS at to finite jj.. The series resummation does not 
diverge, although adding terms causes large changes. The Pade approximant exhibits the divergence but is 
stable under addition of extra terms. 



is a detailed theory of Pade approximants [ |2q ] which needs to be extended to applications in QCD 



where the series coefficients are known only within some statistical errors ||19|]. 

In Figure ^ we show truncated series sums for the QNS at = 0.94Tc. There is no sign of any 
divergence, although successive orders fail to agree with each other as the radius of convergence 
is approached. In the same figure we also show the QNS obtained with Pade resummations of the 
series. These exhibit the divergence identified through series analysis. It is also useful to note that 
the Pade approximants fitted to different number of terms of the series agree with each other except 
when z is significantly larger than the radius of convergence. The Pade analysis indicates a width 
of the critical region which is about A^u/r^ ~ 0.25. 

4. Comparison with experiments 

The QNS are related to fluctuations of conserved quantities in a grand canonical ensemble. It 
may be possible to realize this in an experiment by looking at a part of the fireball produced in a 
heavy-ion collision, provided it thermalizes. Then one way in which grand canonical physics can 
be extracted is by observing particles only in a restricted space-time rapidity range. If this range is 
chosen judiciously, then the remainder of the fireball may act as a heat-bath for the system under 
observation. Then each collision event satisfying the above experimental cuts is one member of a 
grand canonical ensemble of events. 

Event-to-event distributions of conserved quantities then form the observables of interest [^]. 
The cleanest observable is the distribution of total electric charge, Q, since there is very little chance 
of missing a significant fraction of the charge within the acceptance volume V. Baryon number, 
B, and strangeness, S, are also good observables, but since there are uncharged baryons as well as 
long-lived uncharged strange particles, both of which are missed by detectors, the connection to 



these quantities is made at a further remove [29|. Nevertheless, currently the most extensive data 



comes from observations of the net proton number, which is a proxy for the net baryon number. 

It is seen that fluctuations of conserved quantities are Gaussian. The first question is whether 
this Gaussian is entirely (or largely) due to thermal fluctuations. The only way to answer this 
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Figure 9: Comparison of experimental data with lattice predictions. The upper panel shows wi and the 
lower 1112- 

is by going beyond the Gaussian. A systematic way to do this is to change V and check how 
the distribution changes. Gaussian distributions usually arise in experiments through a process 
described by the central limit theorem: with increasing V the higher cumulants of the distribution 



scale down with larger powers of V. The STAR collaboration reported such a measurement [3C] 
using an experimentally determined parameter called the participant number, A'part, as a proxy for 
V. At small A^part all the cumulants, [B"], are comparable, and with increasing A'part the scaling of 
the cumulants is exactly as one should expect — in other words, the microscopic physics encoded 
in the set of [B"] does not change with A'part- 

If the fluctuations are due to thermal processes, then [B"] are related to various NLS com- 
putable in QCD. This is the next step: to check the data against the predictions In order to do 
this one has to take combinations of the cumulants which are free of incidental variables such as 
the uimieasurable V. This is achieved by taking the ratios 

Now, the left hand side of each equation is known from experiment at each ^/S, and the right hand 
side is known from lattice computations if one knows the freezeout values of /i and T. These values 



are parametrized from experimental data in [ |24| ] assuming that the fireball thermalizes. 

on 



On one hand, the series expansions for the NLS are known from the expansion in eq. 3.1 



the other, the ratios of eq. 4.1 have well-determined power behaviour at small z, and poles near 



the critical point. As a result, resummation of the series by Pade approximants is possible. The 
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parameters in the Pade approximants are closely related to the estimates of z*. As a result, the 
lattice artifacts in mi 2,3 are related to those already discussed in the previous section. It turns out 
that both mi and 1112 may have significant finite lattice spacing corrections: mainly a common finite 
multiplicative factor which is also the correction to the estimate to the radius of convergence. For 
m3 lattice spacing corrections are small, except in the vicinity of the critical point. Results for these 



quantities have been given in []33[]. 

The STAR experiment at RHIC has recently published a measurement of comparable quanti- 
ties from runs at three different values of where comparisons with these quantities are presented 



[p2[]. As shown in Figure |9|, it turns out that there is good agreement between the data and the pre- 
dictions. Many questions remain to be answered: on the side of the lattice computations the usual 
questions about flavours, quark masses and lattice spacing, on the experimental side about the 
removal of non-thermal backgrounds and other sources of fluctuations. Nevertheless, this is a sig- 
nificant milestone: the first direct comparison of heavy-ion data with lattice predictions. In future 
such a comparison may even yield a direct measurement of Tc as pointed out in [|^], allowing us 
to set the scale for lattice measurements in an entirely new way. 

5. Assorted topics 

There are various developments at finite chemical potential which cannot be covered fully 
here because of space constraints. However, they are interesting in their own right and have useful 
connections to the physics which is discussed in the previous sections. Here I make a brief mention 
of some of these works. 

In the chiral limit there is a line of second order phase transitions at finite jj., emanating from 
the finite temperature critical point (see Figure [l]). The curvature of this line is an object of interest 
because it sets a scale for the tricritical temperature in the chiral limit. Interesting new work was 



presented on this problem by several groups [34]. 



The phase diagram at imaginary chemical potential is of some interest, since it has to be 
understood if simulations in this region are to be used for analytic continuation to the physically 



interesting case. New results were presented by two groups [|35[]. 

The investigation of correlation functions at finite chemical potential is in its infancy [36|. Any 
new information is interesting at this stage. New work was reported in this meeting [p7|]. 

The strong coupling expansion has been resurrected in this context and improved techniques 
are now being used to investigate the phase diagram. Interesting new results in this direction were 



reported [38]. These may serve to benchmark future simulations using the worm algorithm which 
can be adapted to the strong coupling theory. 

All the results reported till now simulate the grand canonical ensemble. Very little systematic 
effort has gone into simulations in the canonical ensemble with fixed baryon number. One such 



attempt was reported \p%. 
6. Conclusions 

Over the years there has been a great improvement in the understanding of the sign problem at 
finite jJ. in QCD: where it could be tractable and where it is not [Q]. There has been httle progress 
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in directly tackling this problem, although there are interesting developments in the handling of 
other models with sign problems ^ |n|]. However, in the last five years there has been enormous 
progress in lattice computations which can yield information on QCD at finite /i. The essential 
development is the use of analytic continuation, through a Madhava-Maclaurin series expansion in 



/I [13]. There are encouraging tests of this method in performing analytic continuation to finite 



imaginary n where direct simulations are also possible [14,^]. 

The series expansion method has been applied to the extraction of the QCD critical point at 
various lattice spacings [ISj, 19] and with various quark actions and numbers of flavours [|C, 17]. A 



composite figure of the predictions is given in Figure g. It is clear from this figure that the method 
yields results with controlled statistical and systematic errors. Since the results are not strongly 
sensitive to the choice of action, it is also clear that lattice spacing effects are bounded. 

An interesting development in the past year has been the proposal of a set of measurables, eq. 
[4.1| , which allow a direct comparison of experiment and lattice computations. First results show 



very good agreement between data 1 ]32[ ] and prediction [^]. This calls for renewed activity in this 
field and a greater scrutiny of the known systematic uncertainties which need control. 

For communicating their results and then patiently answering my questions I would like to 
thank Gert Aarts, Shailesh Chandrasekharan, Rossella Falcone, Maria-Paola Lombardo, and Chris- 
tian Schmidt. 
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